function  [bdraw,m1,CC1]=drawcon(y,x,prm,prprec,sig,ind,con,convec);
%function  [bdraw,m1,CC1]=drawcon(y,x,prm,prprec,sig,ind,con,convec);
% Draw Coefficients from a Normal distribution when prior is also 
%      Normal( prm, prprec^(-1) )  and prprec is the precision
%      subject to rejection sampling with constraints 
%    
% y      (Txn) 
% x      (TxK)  
% Prm    (Kx1)  Prior Mean 
% prprec (KxK)  Prior Precision
% sig    (1x1)  Variance of the errors
% ind    (Hx1)  indicator of the element numbers to constraint 
% con    (Hx1)  constraints if = -1 less than 
%                              =  1 greater than 
% convec (Hx1)  constraints 
%
% e.g ind=[2 4]   con=[1 -1]  convec=[-10 20]
% restricts the 2nd element to be greater than -10 
%               4th element to be smaller than 20 
%
% bdraw  (1xK)  Draw from the posterior normal as a row vector 
% m1     (1xK)  Posterior mean 
% CC1    (KxK)  Cholelsky of the posterior variance 
%               V1= CC1'*CC1  
[ty,ny]=size(y); 
[tx,nx]=size(x); 
nk=length(prm);
if tx~=ty
    error('Not equal number of observations in drawcoef'); 
end; 

V1=( prprec+ (1/sig)*(x'*x)  ) \ eye(nk) ;   % Posterior Variance  
m1=V1*(prprec*prm + (1/sig)*(x'*y) );        % Posterior Mean  
CC1=chol(V1); 
cons=con.*convec;

acc=0;  
rcoun=1; 

while acc <1
    bdraw = m1+ CC1*randn(nk,1) ;
    if ( con .* bdraw(ind)  ) > cons; 
        acc=1;  
    end 
    rcoun=rcoun+1; 
    if rem(rcoun,500) == 0; 
        warning('500 draws have not yet satisifed the constraint in drawNcon.m') 
    end 
end 
bdraw=bdraw'; 
m1=m1';
rcoun=rcoun-1; 